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We present an improved upper bound for the ground state energy of lattice fermion models with 
sign problem. The bound can be computed by numerical simulation of a recently proposed family 
of deformed Hamiltonians with no sign problem. For one dimensional models, we expect the bound 
to be particularly effective and practical extrapolation procedures are discussed. In particular, in 
a model of spinless interacting fermions and in the Hubbard model at various filling and Coulomb 
repulsion we show how such techniques can estimate ground state energies and correlation function 
with great accuracy. 



In the study of strongly interacting quantum systems, 
the determination of ground state properties is a basic 
step in the analysis of many physically relevant prob- 
lems. To accomplish this task, several powerful numeri- 
cal methods are available and big efforts have been spent 
in their development during the last years However, 
as is well known, in fermion models, numerical simula- 
tions face serious obstructions that dramatically reduce 
their performance In this Brief Report, we propose 
an improved upper bound for the ground state energy 
and an effective strategy to circumvent this difficulty in 
a certain class of one dimensional models. 

For moderate lattice size, Exact Diagonalization meth- 
ods are possible ||^ with no additional problems in the 
fermionic case. Typically, these are calculations based on 
the Lanczos algorithm and the ground state of systems 
with up to several thousands states can be determined. 
For larger state space, other methods must be consid- 
ered. In particular, the most flexible tool is Monte Carlo 
simulation that evaluates quantum expectation values by 
a careful stochastic sampling of the configuration space. 

Often, when the model under study is fermionic, the 
measure to be sampled is non-positive. This fact can 
be related to the intrinsic anticommuting nature of the 
Fermi creation annihilation operators, but is more gen- 
eral and appears also in bosonic contexts, like quantum 
spin models with non trivial exchange terms Q . This sit- 
uation, usually described as the "sign-problem" , is quite 
serious and standard algorithms simply fail. 

Several proposals are available to overcome or at least 
reduce the sign problem in specific problems, in partic- 
ular realistic models in more than one dimension. One 
of the most powerful techniques is the Constrained Path 
Monte Carlo Method where the sign problem is elimi- 
nated by introducing a physically motivated guiding wave 
function that approximates the nodal structure of the 
exact ground state. A numerical analysis is then pos- 
sible and expectation values can be computed. The 
approximation is partially controlled and, for instance. 



the estimates of the ground state energy are known to 
provide rigorous upper bounds. Other alternatives are 
the Fixed-Node approximation, the Projection Quantum 
Monte Carlo, path integral representations or the Aux- 
iliary Field Monte Carlo Of course, trial wave func- 
tions are useful also in cases with no sign problem, but 
their purpose to accelerate Monte Carlo simulations and 
the final results are in principle independent on the choice 
of the guiding function. 

Our first goal will be that of providing an improved 
general upper bound for the energy of fermion lattice 
models with sign problem; the bound is computable by 
numerical simulations. 

As a second aim, we ask whether one dimensional prob- 
lems allow for any special simplification and practical 
treatment of the sign problem in the spirit of 0. In 
particular we look for methods that do not rely on any 
biased approximation. A hint to a positive answer comes 
from the remark that in one dimension, the sign prob- 
lem appears to be somewhat artificial. For instance, we 
can consider a model where the standard fermion hop- 
ping term is the unique source of alternating signs in the 
off diagonal matrix elements of the Hamiltonian. In this 
case, the sign problem disappears as soon as open bound- 
ary conditions are used |^ . Its manifestations are there- 
fore expected to be completely negligible in observable 
quantities that admit a thermodynamical limit, indepen- 
dent on the boundary conditions. Possible examples are 
the ground state energy per site and finite correlation 
lengths. For this reason, we are led to believe that in such 
cases, the finite size effect of the sign problem should be 
under control without much effort. 

To begin our analysis, let us consider a d dimensional 
hypercubic lattice with V — L1L2 ■ ■ ■ Ld sites and peri- 
odic boundary conditions. Let us assign to each site i a 
pair of spinless fermion creation annihilation operators Ci 
and C; with canonical algebra 

= {ci,cj} = 0, {4,c]} = 0. (1) 
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Let us denote a configuration by n, the set of occupa- 
tion numbers at all the lattice sites. We work in the 
occupation number basis with states {|n)} such that 
c[ci|n) ~ 7ii|n) and assume a given ordering of the sites 
to fix signs. We want to study Hamiltonians of the form 



H — Hkii- 



Hi, 



where 



Hkin — 



(iJ> 



cjci), 



(2) 



(3) 



((i, j) denotes a pair of connected sites, for instance near- 
est neighbours). The interaction term is taken to be di- 
agonal in the { | n) } basis 



Hi, 



y(n) 



(4) 



with V being a real function. The fermion number 
TV = rii is conserved and the analysis of the spectrum 
can be carried on in each sector with definite N. 

In the one dimensional case, V = Li = L, we can 
assume the following canonical site ordering 



0), 



(5) 



where |0) is the empty state. The paired sites are {i, 
modulo periodic boundary conditions. The off diagonal 
matrix elements of H are then negative and Monte Carlo 
simulations do not have sign problems, except when the 
fermion number N is even; in this case, there exist ma- 
trix elements with the wrong sign, precisely those con- 
necting pairs of states that differ by a hopping through 
the boundary. 

As shown in 1^ , it is possible to introduce a new Hamil- 
tonian, _ff°*^, with no sign problem. Its ground state en- 
ergy Eq^ can be computed by numerical methods and 
provides a rigorous upper bound for the ground state en- 
ergy Eq of H. In this result has been extended by 
introducing a family of Hamiltonians {ff (7)} depending 
on a real free parameter 7. Let us recall the definition 
of ^(7) = {Hijij)} in terms of iJ = {Hij}- We assume 
H = DH^^^D^^ with H^°'^ hermitean (for instance of 
the form Eq. (||)) and D diagonal with Da 7^ 0. D is 
a trial wave function useful to accelerate convergence to 
the ground state. The new Hamiltonian i?(7) is given by 




H, 
H, 



>0,k^i 



H 



ik 



< i 
>0 z^j 

i = 3 

(6) 



If we call Ef) (7) the ground state energy of (7) , then it is 
possible to show that: (i) Eo{—l) = Eq, (ii) Eo{j) > Eq 
for any 7 > — 1, (iii) H{'-f) has no sign problem if 7 > 0. 
In the Monte Carlo approach adopted in (nOl, it is also 



emphasized that a strictly positive 7 must be used and 
that an extrapolation to 7 = is required for best accu- 
racy. 

We now show that the function Eq{"/) is concave in 7. 
Namely, for < a < 1 and any two real 71 , 72 , we have 

So(a7i + (1 - a)72) > a£^o(7i) + (1 - a)£'o(72)- (7) 

A general simple proof can be given by means of the vari- 
ational characterization of Eq. The dependence of H{'j) 
on the parameter 7 is linear and we can write 



H{"f) = DH^^\-i)D-^ = D{Hi +-fH2)D-\ 



(8) 



with hermitean, 7 independent Hi and i?2- Then, if S is 
the set of normalized stated, we have 



£^o(a7i + (1 - ")72) 



min 



(9) 



> amin(s|i7(")(7i)|s) + (1 
a-Bo(7i) + (1 - a)So(72). 



a) min(s 



\H^'\l2M 



Concavity means that the incremental ratio is a decreas- 
ing function. Computing the incremental ratio on the 
intervals [—1,0] and [71,72] with arbitrary < 71 < 72 
we find 



Eq < Eo{0) 



£^0(72) - £^0(71) 



72 - 71 



(10) 



If we furthermore assume i?o(7) to be differentiable at 
7 = and take the limit 72,71 ^ we can replace the 
the upper bound Eq < Eq(0) by the the improved one 



Eq < Eq{0) - E'q{0^ 



(11) 



that expresses the geometrical fact that a concave curve 
lies below its tangent at any point. We remark that Eq{j) 
can have finite jumps. If they can be excluded a priori, 
than the above proof follows also by using second order 
perturbation theory. 

The concavity constraint suggests that the function 
i?o(7) could be particularly well behaved. For this rea- 
son we also try to obtain an extrapolation of its value 
at 7 = —1 by assuming it smooth enough. To be fair, 
this possibility must be considered just a practical recipe. 
Nonetheless, we shall give numerical arguments to sup- 
port its robustness. 

We begin with a simple model of interacting spinless 
fcrmions. The Hamiltonian is 



H 



(12) 



V 



H+l 



For lattice size L — 8, 12, 16 and 20, at half filling and for 
several interaction couplings V, we determine by Lanc- 
zos diagonalization the exact ground state energy Eq, the 
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bound Eg (0) , the derivative rj = E'q (0) and the improved 
bound Eq"^^^ = Eo{0) — t]. Moreover, we also attempt 
an extrapolation to 7 = — 1 starting from the values of 
i?o(7) in the range 7 £ (0, 1). In principle, any such ex- 
trapolation can be very dangerous. On the other hand, 
the extrapolated Eq"^\—1) obtained by a polynomial fit 
of degree m displays a flat behaviour with small oscilla- 
tions unless m is too large when the extrapolation process 
breaks down. The typical degree of the best polynomial 
is 8. To accelerate the convergence and to determine the 
best estimate, we also used Aitkin's algorithm that 
improves the converging sequence {zn} replacing it with 
{z^} defined by 

' — _ (AZn) ^ ^ (TO^ 

Z^ — Zn A 2 ' — ^n; \^^} 

(where A^z„ = A(Az„)). We call E^f^*'^ the resulting 
prediction. In Fig. 1, we show the results for Eq, Eo{Q), 
j^j^mp) in terms of their relative percentual ac- 

curacy defined as 100 • \5E/Eo\ where 5E is the error in 
the ground state energy estimate. The improved bound 
is significantly better than i?o(0). Both -Eo(O) and iJ^*™^' 
converge to the exact value as the system size L — » 00 
as expected from our initial discussion on the asymptotic 
irrelevance of the boundary effects. The extrapolated 
bound e'^^^^ is quite precise and around the permille 
level. 

A similar analysis can be done for the measurement of 
observables. Unfortunately, for these, we cannot derive a 
simple bound like the one for the energy and in fact ran- 
dom hermitean operators produce easily wild behaviour 
for 7 G (—1,0). On the other hand, for operators as- 
sociated to physically meaningful quantities independent 
on the boundary conditions we can expect a mild depen- 
dence on 7. As a typical non local observable, we consider 
the integrated staggered correlation function defined by 

-, L L/2 
1=1 k=l 

In Fig. 2, we show how the exact values are very accu- 
rately reproduced by the extrapolated values. The rela- 
tive accuracy is always well below the percent level. In 
Fig. 3, we show the behaviour of the functions i?o(7) and 
5(7) aXV = l for the system with L — IQ together with 
a 8-th degree polynomial fit to emphasize their smooth- 
ness. 

For large systems, the upper bound must be computed 
by numerically extrapolating -Eo(O) and Eq{0). This 
can be done by a straightforward Monte Carlo simu- 
lation at variable 7. To explore the feasibility of this 
proposal, we perform such a study for free fermions 
{V = 0) on a lattice with L = 40 sites at half fill- 
ing. Here, the exact ground state energy is known, 



Eo/2t = — cot(7r/L) —12.706 and can be used as a 
check. We use the Green Function Monte Carlo with 
Stochastic Reconfiguration ||l^.|l^. We run simulations 
with variable number of walkers in order to extrapolate 
the infinite population size limit. In Fig. 4 we show the 
estimated value of Eq (7) at several positive 7 values to- 
gether with a simple parabolic fit. The estimated value 
of the bound is EoiO) - E'„{0) = -12.61(1), about 1% off 
the exact value. 

Similar results are obtained by studying the one di- 
mensional Hubbard model. Denoting by tj i the two 
spin degrees of freedom, the Hamiltonian reads 

H = -t (4<TCi+i,<T + cl+i,aCt,a) + UY^ ^^i,T"^i• 

'i>o"=t,X i 

(15) 

Again, we determine for several lattice size, filling frac- 
tion and coupling U the four quantities Eq, Eo{0), E^^^ 
and £^Q^'. The results are shown in Fig. 5 where we show 
the relative percentual accuracy of Eq'^^*\ The quality 
of the results is not good as in the spinless model, but 
the errors are again small, of a few percents. In fact, a 
scaling analysis shows convergence to the exact values in 
the large volume limit. 

To summarize, the family of Hamiltonian with no sign 
problem proposed in |ic| ] makes possible the derivation of 
a size consistent bound for the ground state energy that 
improves the one at 7 = 0. Moreover, much information 
can be reconstructed for the original Hamiltonian. The 
accuracy of our analysis on small systems is certainly be- 
yond a practical implementation, but suggests that also 
for more complicated systems, not allowing a direct anal- 
ysis, the extrapolated values can provide useful numer- 
ical hints. Preliminary results on the two dimensional 
Hubbard model are encouraging and will be presented 
elsewhere JT^ . 
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FIG. 2. Spinless interacting fermions. Relative percentual 
accuracy of the extrapolated estimate of the integrated stag- 
gered correlation function S. 
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FIG. 1. Spinless interacting fermions. Relative percentual 
accuracy of the ground state energy bounds and extrapolated 
estimate. The circles with very small error are the best ex- 
trapolated bound Eq'^*^ at the various lattice sizes considered 
L = 8, 12, 16, 20 (we use the same symbol since the aim is to 
emphasize the small maximum deviation). 
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FIG. 3. Spinless interacting fermions. Smoothness of the 
function £0(7), 'S'(7) and accuracy of a 8-th degree polyno- 
mial fit for the model at \^ = 1 with L = 16 sites. 
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FIG. 4. Spinless interacting fermions. GFMC calculation 
of the effective energy i?o(7) in the free case V = 0, t = 0.5, 
on a lattice with L = 40 sites at half filling. Data are obtained 
by the GFMC with Stochastic Reconfiguration in the Poisson 
limit (continuous time). The size of the walker ensemble has 
been varied from 300 up to 2000 and the extrapolation to 
the infinite population size has been taken. Statistics corre- 
sponds to 10^ reconfigurations separated by an evolution time 
T = 0.5. The dashed line is a parabolic fit. 
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FIG. 5. Hubbard model. Relative percentual accuracy of 
the extrapolated estimate of the ground state energy at dif- 
ferent sizes, filling and Coulomb coupling. 
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